Análise de dados

Com ajuda do tidytuesday (dados de hospedagem em hotel)

Ianní Muliterno
Geovana Lopes

Ministrante

  • Ianní Muliterno
  • Graduado em Estatística na Universidade Federal de Pernambuco
  • Cientista de dados na Unilever
  • Co-organizador da R-Ladies São Paulo, uma comunidade que tem como objetivo promover a diversidade de gênero na comunidade da linguagem R.
  • Escritor
  • Geovana
  • Graduada em Matemática (licenciatura). Atualmente cursando Estatística (bacharel) no IME-USP e Pós-Graduação em Ciência de Dados no IFSP Campinas;
  • Cientista de Dados na Almap BBDO;
  • Co-organizadora da R-Ladies São Paulo.

Pré-requisitos

  • R e RStudio instalados no seu computador:

  • Links para instalação:

Análise de dados

A análise de dados tem por objetivo transformar informação em insights e soluções, não necessariamente envolve modelagem, mas gosto muito desta visão do ciclo da análise de dados:

:::footer Créditos da imagem Allison Horst. :::

A nossa base

A base de hoje se chama hotels, compartilhada pelo Tidytuesday em 2020. Tidytuesday é um projeto que lança semanalmente bases de dados reais. Organizado pela comunidade R4DS.

Primeiro passo - Importar

Quando a base de dados está disponível em algum repositório público do github, podemos carregar assim:

#importando pacotes
library(tidyverse)
library(skimr)
library(GGally)
library(tidymodels)
library(recipes)
library(workflows)
library(RColorBrewer)
library(ggcorrplot)


hotels <- readr::read_csv("https://raw.githubusercontent.com/rfordatascience/tidytuesday/master/data/2020/2020-02-11/hotels.csv")

Um problema comum que podemos encontrar é ter de carregar uma base de dados local muito grande, neste caso, uma boa saída seria utilizar data.table

#base <- data.table::fread("seu_path")

Segundo passo - Arrumar (Limpeza e pré-processamento)

Primeiramente, precisamos conhecer a base,

  • Quais colunas temos disponíveis?
  • Quais os tipos das colunas?
  • Existe algum problema de coleta? . . . Aqui temos o jeito prático e o jeito amigável. O prático serve para o analista, o amigável para seu consumidor final.

Arrumar - Prático vs Amigável

Para a versão prática, vamos utilizar a função skim do pacote skimr, para obtermos um overview.

Prático

Prático

Prático

Amigável

hotels |> 
  select_if(is.numeric) |> 
  pivot_longer(1:18) |> 
ggplot(aes(x=value)) +
  geom_histogram() +
  labs(title="Histograma") +
  facet_wrap(~name,scales = "free")

Amigável

Arrumar - O que poderia ter acontecido

Lidar com dados vazios, inputs errados e muitas NAs, geralmente levam a:

  • Contato com o fornecedor dos dados (pode ser interno, caso exista um time de administração de banco de dados, ou data engineering), para entender a possibilidade de correção.

  • Remoção das colunas ou linhas afetadas (preencher dados faltantes não é algo muito utilizado na prática, por ser arriscado)

Terceiro passo - Transformar

Considerar apenas as estadias que de fato aconteceram pois nosso objetivo na análise envolve os hóspedes com filhos.

hotel_stays <- hotels |> 
  filter(is_canceled == 0) |> 
  mutate(
    children = case_when(
      children + babies > 0 ~ "children",
      TRUE ~ "none"
    ),
    required_car_parking_spaces = case_when(
      required_car_parking_spaces > 0 ~ "parking",
      TRUE ~ "none"
    )
  )  |> 
  select(-is_canceled, -reservation_status)

Quarto passo - Visualizar

Chegou a hora de fazer perguntas para os dados.

hotel_stays |> 
  mutate(arrival_date_month = factor(arrival_date_month,
                                     levels = month.name
  )) |> 
  count(hotel, arrival_date_month, children) |> 
  group_by(hotel, children) |> 
  mutate(proportion = n / sum(n)) |> 
  ggplot(aes(arrival_date_month, proportion, fill = children)) +
  geom_col(position = "dodge") +
  scale_y_continuous(labels = scales::percent_format()) +
  facet_wrap(~hotel, nrow = 2) +
  labs(
    x = NULL,
    y = "Proporção de estadia ao longo do ano",
    fill = NULL
  )

Visualizar

Visualizar

hotel_stays |> 
  count(hotel, required_car_parking_spaces, children) |> 
  group_by(hotel, children) |> 
  mutate(proportion = n / sum(n)) |> 
  ggplot(aes(required_car_parking_spaces, proportion, fill = children)) +
  geom_col(position = "dodge") +
  scale_y_continuous(labels = scales::percent_format()) +
  facet_wrap(~hotel, nrow = 2) +
  labs(
    x = NULL,
    y = "Proporção de estadia",
    fill = NULL
  )+
  ggtitle("Proporção de acordo com solicitação de vagas e presença de crianças")

Visualizar

Visualizar

Características de hospedagem

df_processed <- hotels |> 
  filter(is_canceled == 0) |> 
  mutate(
    parents = if_else(children > 0 | babies > 0, 1, 0),
    parents_detailed = case_when(children > 0 ~ "pais",
                                 babies > 0 ~ "pais_jovens",
                                 TRUE ~ "sem_filhos"),
    meal_num = if_else(meal %in% c("HB", "FB"), 1, 0)
  ) |> 
  dplyr::select(parents_detailed,# not_canceled,
                stays_in_weekend_nights, 
                meal_num, total_of_special_requests) |> 
  pivot_longer(cols = stays_in_weekend_nights:total_of_special_requests,
               names_to = "variable",
               values_to = "value") |> 
  group_by(parents_detailed, variable) |>
  summarize(
    yes = sum(value > 0) / n(),
    no = sum(value == 0) / n()
  ) |> 
  pivot_longer(cols = c(yes, no),
               names_to = "group",
               values_to = "value")|>
  filter(group == "yes")

Visualizar

Características de hospedagem

# Paletas de cores simplificadas
color_fill <- c("pais_jovens" = "#66c2a5", "pais" = "#fc8d62", "sem_filhos" = "#8da0cb")

  ggplot(df_processed , aes(x = variable, y = value, fill = parents_detailed)) +
  geom_bar(stat = "identity", position = "dodge") +
  scale_fill_manual(values = color_fill) +
  scale_alpha_manual(0.6) +
  labs(title = "Comparação de característica de hospedagem",
       x = "Características",
       y = "%",
       fill = "tem filhos?",
       alpha = "Type of Parent") +
  theme_minimal()

Visualizar

Características de hospedagem

Antes de modelar…

primeiramente, quais são as variáveis mais relacionadas com a variável resposta?

Ggpairs nos ajuda a analisar de modo geral o que tem correlação com a presença de crianças

hotel_stays |>
  select(
    children, adr,
    required_car_parking_spaces,
    total_of_special_requests
  ) |>
  ggpairs(mapping = aes(color = children))

Visualizar

pausa para apreciar o poder dos histogramas+boxplot

Mas…

  • Correlação Clássica (Pearson): Esta é utilizada para avaliar a relação linear entre duas variáveis contínuas ou ordinais. Assim, é apropriada para variáveis numéricas.

  • Correlação de Spearman: Esta é uma medida de correlação de postos e não pressupõe uma relação linear. É frequentemente utilizada para variáveis ordinais, mas também pode ser usada com variáveis numéricas quando a relação entre elas não é estritamente linear.

  • Qui-quadrado: Serve para ter noção se variáveis categóricas estão transmitindo informações semelhantes, algo menos usual que eu particularmente não aplico no dia a dia se preciso construir modelo.

Outras questões

  • Também poderíamos ter de lidar com multicolinearidade
  • Se a base fosse muito grande poderíamos precisar rodar uma seleção de variáveis

E aí?

Multicolinearidade X Seleção de variáveis

  • Rodar uma seleção de variáveis (por exemplo, stepwise), não garante que a multicolinearidade será tratada, então ambos devem ser feitos em caso de um modelo muito grande
  • Dito isto, qual o tamanho do problema que é, ter multicolinearidade?
  • como evitar ? (uma forma nós ja vimos)

VIF

Avalia se existe relação múltipla entre variáveis X’s. A pergunta que ela responde é:

  • o quanto o \(R^2\) da variavel \(X_1\) assumindo o papel de Y com todas as outras X´s?

O pacote olsrr oferece a função ols_vif_tol(modelo)

Quinto passo - Modelar

Selecionar as colunas mais importantes para a variável resposta e transformar as categóricas em fator (boa prática, pois alguns modelos, como knn podem demandar isso).

hotels_df <- hotel_stays |>
  filter(adr > 0) |> 
  select(
    children, hotel, arrival_date_month, meal, adr, adults,
    required_car_parking_spaces, total_of_special_requests,
    stays_in_week_nights, stays_in_weekend_nights
  ) |>
mutate_if(is.character, factor)

Modelar

set.seed(2023)

#particionar os dados em treino e teste
hotel_split <- initial_split(hotels_df)

hotel_train <- training(hotel_split)
hotel_test <- testing(hotel_split)

# Cria uma "receita" para o pré-processamento dos dados. 
# Estamos tentando prever a variável 'children' com base em todas as outras variáveis (denotadas por '.') no conjunto de dados 'hotel_train'.
hotel_rec <- recipe(children ~ ., data = hotel_train) |>
  #  método de subamostragem para balancear a variável de resposta 'children'. 
  themis::step_downsample(children) |>
  # Converte todas as variáveis nominais (categóricas) em variáveis dummy (ou binárias), exceto a variável de resposta. 
  # Por exemplo, uma variável categórica com níveis 'A', 'B', e 'C' será convertida em três novas variáveis binárias.
  step_dummy(all_nominal(), -all_outcomes()) |>
  # Remove qualquer variável numérica que tenha um único valor (variância zero) em todas as observaçõesjá que tais variáveis não adicionam informações úteis ao modelo.
  step_zv(all_numeric()) |>
  # Normaliza todas as variáveis numéricas para que tenham média 0 e desvio padrão 1. 
  # Isso ajuda em muitos algoritmos de aprendizado de máquina que são sensíveis à escala das variáveis.
  step_normalize(all_numeric()) |>
  # Prepara a receita aplicando todas as etapas definidas acima.
  prep()

hotel_rec

É possível usar recipes apenas para pré-processamento e então partir para as funções de modelo que estamos acostumados. (exemplo)

Modelar - GLM

# Especificar o modelo GLM
glm_spec <- logistic_reg() |>  
  set_engine("glm") |>  
  set_mode("classification")

# Ajuste o modelo
glm_fit <- glm_spec |>  
  fit(children ~ ., data = juice(hotel_rec))
glm_fit
parsnip model object


Call:  stats::glm(formula = children ~ ., family = stats::binomial, 
    data = data)

Coefficients:
                        (Intercept)                                  adr  
                           -0.08878                             -1.32831  
                             adults            total_of_special_requests  
                            0.10982                             -0.45245  
               stays_in_week_nights              stays_in_weekend_nights  
                           -0.01226                             -0.07458  
                 hotel_Resort.Hotel            arrival_date_month_August  
                           -0.12503                             -0.01568  
        arrival_date_month_December          arrival_date_month_February  
                           -0.15141                             -0.18606  
         arrival_date_month_January              arrival_date_month_July  
                           -0.10229                             -0.04909  
            arrival_date_month_June             arrival_date_month_March  
                            0.13670                             -0.03547  
             arrival_date_month_May          arrival_date_month_November  
                            0.11576                             -0.01987  
         arrival_date_month_October         arrival_date_month_September  
                            0.02699                              0.21008  
                            meal_FB                              meal_HB  
                            0.04334                              0.17046  
                            meal_SC                       meal_Undefined  
                            0.25721                              0.08655  
required_car_parking_spaces_parking  
                           -0.12133  

Degrees of Freedom: 8919 Total (i.e. Null);  8897 Residual
Null Deviance:      12370 
Residual Deviance: 9496     AIC: 9542

Modelar - GLM

validation_splits <- mc_cv(juice(hotel_rec), prop = 0.9, strata = children)
# Definir o fluxo de trabalho para o GLM
glm_wf <- workflow() |>  
  add_model(glm_spec) |>  
  add_formula(children ~ .)



glm_res <- fit_resamples(
  glm_wf,
  validation_splits,
  control = control_resamples(save_pred = TRUE)
)

glm_res |>
  collect_metrics()
# A tibble: 2 × 6
  .metric  .estimator  mean     n std_err .config             
  <chr>    <chr>      <dbl> <int>   <dbl> <chr>               
1 accuracy binary     0.748    25 0.00211 Preprocessor1_Model1
2 roc_auc  binary     0.817    25 0.00238 Preprocessor1_Model1

Modelar - KNN

knn_spec <- nearest_neighbor() |>
  set_engine("kknn") |>
  set_mode("classification")

knn_fit <- knn_spec |>
  fit(children ~ ., data = juice(hotel_rec))

knn_fit
parsnip model object


Call:
kknn::train.kknn(formula = children ~ ., data = data, ks = min_rows(5,     data, 5))

Type of response variable: nominal
Minimal misclassification: 0.2543722
Best kernel: optimal
Best k: 5

Modelar - KNN

knn_wf <- workflow() |> 
  add_model(knn_spec) |> 
  add_formula(children ~.)

knn_res <- fit_resamples(
  knn_wf,
  validation_splits,
  control = control_resamples(save_pred = TRUE)
)

knn_res |>
  collect_metrics()
# A tibble: 2 × 6
  .metric  .estimator  mean     n std_err .config             
  <chr>    <chr>      <dbl> <int>   <dbl> <chr>               
1 accuracy binary     0.742    25 0.00262 Preprocessor1_Model1
2 roc_auc  binary     0.808    25 0.00297 Preprocessor1_Model1

Modelar - Árvore de decisão

tree_spec <- decision_tree() |>
  set_engine("rpart") |>
  set_mode("classification")

tree_fit <- tree_spec |>
  fit(children ~ ., data = juice(hotel_rec))

tree_fit
parsnip model object

n= 8920 

node), split, n, loss, yval, (yprob)
      * denotes terminal node

 1) root 8920 4460 children (0.5000000 0.5000000)  
   2) adr>=-0.06730205 3982 1041 children (0.7385736 0.2614264) *
   3) adr< -0.06730205 4938 1519 none (0.3076144 0.6923856)  
     6) total_of_special_requests>=0.6624924 890  406 children (0.5438202 0.4561798) *
     7) total_of_special_requests< 0.6624924 4048 1035 none (0.2556818 0.7443182)  
      14) adults< -2.89155 63    0 children (1.0000000 0.0000000) *
      15) adults>=-2.89155 3985  972 none (0.2439147 0.7560853) *

Modelar - Árvore de decisão

tree_wf <- workflow() |> 
  add_model(tree_spec) |> 
  add_formula(children ~.)


tree_res <- fit_resamples(
  tree_wf,
  validation_splits,
  control = control_resamples(save_pred = TRUE)
)

tree_res |>
  collect_metrics()
# A tibble: 2 × 6
  .metric  .estimator  mean     n std_err .config             
  <chr>    <chr>      <dbl> <int>   <dbl> <chr>               
1 accuracy binary     0.733    25 0.00283 Preprocessor1_Model1
2 roc_auc  binary     0.759    25 0.00402 Preprocessor1_Model1

Comparação de modelos

knn_res |> 
  unnest(.predictions) |> 
  mutate(model = "kknn") |> 
  bind_rows(tree_res |> 
              unnest(.predictions) |> 
              mutate(model = "rpart")) |> 
  bind_rows(glm_res |> 
              unnest(.predictions) |> 
              mutate(model = "glm")) |> 
  group_by(model) |> 
  roc_curve(children, .pred_children) |> 
  ggplot(aes(x = 1 - specificity, y = sensitivity, color = model)) +
  geom_line(size = 1.5) +
  geom_abline(
    lty = 2, alpha = 0.5,
    color = "gray50",
    size = 1.2
  )

Comparação de modelos

Melhor modelo + base de teste

Outra coisa incomum no meio acadêmico que é comum nos negócios é separar em treino e teste para trazer mais segurança sobre a performance do modelo

test_proc <- bake(hotel_rec, new_data = hotel_test)

glm_fit |> 
  predict(new_data = test_proc, type = "prob")  |> 
  mutate(truth = hotel_test$children) |> 
  roc_auc(truth, .pred_children)
# A tibble: 1 × 3
  .metric .estimator .estimate
  <chr>   <chr>          <dbl>
1 roc_auc binary         0.811

Agradecimentos

R-Ladies theme for Quarto Presentations. Code available on GitHub.